Two-Flux Green’s Function Analysis 
for Transient Spectral Radiation 
in a Composite 

R. Siegel 


Reprinted from 

Journal el Thrnnopliysies ami Heat Transfer 

Volume 10, Number 4, Pages 681-688 


A publication of the 

American Institute of Aeronautics and Astronautics, Inc. 
1801 Alexander Bell Drive, Suite 500 
Reston, VA 22091 




Journal of Thermophysics and Heat Transfer 
Vol. 10, No. 4, October- December 1996 


Two-Flux Green’s Function Analysis for Transient Spectral 

Radiation in a Composite 


Robert Siegel* 

NASA Lewis Research Center, Cleveland, Ohio 44135 


An analysis is developed for obtaining transient temperatures in a two-layer semitransparent composite 
with spectrally dependent properties. Each external boundary of the composite is subjected to radiation 
and convection. The two-flux radiative transfer equations are solved by deriving a Green’s function. This 
yields the local radiative heat source needed to numerically solve the transient energy equation. An 
advantage of the two-flux method is that isotropic scattering is included without added complexity. The 
layer refractive indices are larger than one. This produces internal reflections at the boundaries and the 
internal interface; the reflections are assumed diffuse. Spectral results using the Green’s function method 
are verified by comparing with numerical solutions using the exact radiative transfer equations. Transient 
temperature distributions are given to illustrate the effect of radiative heating on one side of a composite 
with external convective cooling. The protection of a material from incident radiation is illustrated by 
adding scattering to the layer adjacent to the radiative source. 
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Nomenclature 

= spectral absorption coefficient of each layer, trT 1 
= coefficients in Green’s functions 
= specific heat of each layer in the composite, 
J/kg-K 

= speed of electromagnetic propagation in 
vacuum, m/s 

= heat capacity ratio of layers, p 2 c 2 /p t c 1 
= total thickness of composite, D x + D 2 , m 
= thickness of each layer in composite, m 
= coefficients in homogeneous solution for 
each layer 

= blackbody spectral energy in vacuum divided 
by crT ? 

= quantity in each layer, k Jv (D/DJ) [3(1 — fX,v)] 1/2 
= flux quantity 2 (q? + q7), W/m 2 ; 

G = G/crT 4 

= Green’s function for G, in each layer 
= dimensionless parameters, hJcrT 3 and 
hJo-T] 

- convective heat transfer coefficients at x — 0 
and D, W/m 2 K 

= extinction coefficient of each layer, aj + cr sj , 
m 

=s thermal conductivity of each layer, W/m • K 
= conduction -radiation parameter of each 
layer based on D, kj!4oT]D 
= refractive index of each layer 
= quantities P 2v — 2R 7 !3K iv D, 

P Av = 2R 4 /3K 1 v D, P 6v = 2RJ3K 2v D 
= radiative flux in the x direction, W/m 2 ; 
q r = qJcrT 4 

= external radiation fluxes crT* x and crT % 
incident at x — 0 and D, W/m 2 
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= dimensionless radiation fluxes, q rl /crT* and 
qn/a-H 

= radiative fluxes in positive and negative x 
directions, W/m 2 

= reflection coefficients, = (1 — p j/( 1 — p b ), 
R 2 = (1 + ft,)/( 1 - PtX /?3 = 0 “ PcV 

(1 - Pd ), i? 4 = (Pc + Prf)/(1 ” P«\ Rs = 

(1 - Py)/( 1 — Pe), R(, ~ (1 + p7)/{\ ~ Pe ) 

= absolute temperature, K; t = TITi 
= gas temperatures for convection at x = 0 and 
D, K 

= reference temperature, or uniform initial 
temperature, K 

= temperatures of blackbody radiative 
surroundings at x = 0 and D, K; 

T sl = {q rl /o-y'\ T s2 = (q r2 /cr) U4 
= dimensionless gas temperatures, T xl /T„ T s2 IT t 
= coordinate in direction across layer, m; 

X = x/D 

~ dimensionless thickness of first layer, DJD 
= time, s 

= optical thickness of each layer, (a jv + cr v JD y 
= frequency of spectral radiation, s~ l 
= cutoff frequency in each layer, s~‘ 

= location within composite divided by D 
— reflectivities at interfaces. Fig. 1 
= density of each layer, kg/m 3 
= Stefan -Boltzmann constant, W/m 2 K 4 
= spectral scattering coefficient of each layer, 
m _1 

= dimensionless time, (4cr7’ 3 /p 1 c J D)0 
= spectral scattering albedo of each layer, 

CT sj J) 


Subscripts 

h — homogeneous solution 

j = layer index, j = 1 in first layer, j - 2 in 

second layer 


Introduction 

D ETAILED transient temperatures are needed to examine 
heat transfer and thermal stress behavior for materials at 
high temperature that are partially transparent to thermal ra- 
diation. Applications include ceramic components, thermal 
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protection coatings, tempering of glass windows, high-temper- 
ature insulation systems, and porous ceramic burners. For 
semitransparent materials where thermal radiation affects in- 
ternal temperatures, transients have been studied much less 
than steady state. 

For transient solutions, numerical procedures such as finite 
differences and finite elements have been used to solve the 
transient energy equation coupled with the radiative transfer 
relations. Some of the literature is reviewed in Refs. 1—7, 
and so a detailed review of transient solutions is not repeated 
here. A few investigations of particular interest for the present 
study provide analyses of composites with two or multiple 
layers. 8- " In Ref. 8, the radiative transfer equations for ab- 
sorbing layers with diffuse boundaries are solved numerically 
in combination with the transient energy equation; scattering 
is not included. A finite difference method is used in the spatial 
direction, and a Runge— Kutta method is used to move forward 
in time. A multilayered region with diffuse boundaries is con- 
sidered in Ref. 9; scattering is not included. A finite difference 
method is used with a variable grid size. A two-layer com- 
posite is analyzed in Ref. 10 subjected to external radiation at 
one of its boundaries. Isotropic scattering is included, and the 
layers are gray. A Galerkin method is used to solve the radi- 
ative transfer relations, coupled with a finite difference solution 
of the transient energy equation. 

The present analysis continues the development in Ref. 1 1 
where a two-flux Green’s function method was developed for 
a two-layer gray composite; spectral properties are now in- 
cluded. The two-flux approximation 12 is used for the spectral 
radiative transfer equations. An advantage of the two-flux 
method is that isotropic scattering is included without addi- 
tional effort, as compared with the need to solve the scattering 
integral equation at each time when using the exact transfer 
relations. For the general boundary conditions of external con- 
vection and radiation for a layer having diffuse interfaces, it 
was shown in Refs. 13 and 14 that the two-flux method pro- 
vides accurate temperature distributions and heat fluxes. This 
was for steady state in single and composite layers with scat- 
tering included. For transient behavior, the two-flux method 
was shown to provide accurate results for a single layer, 3 and 
for a two-layer gray composite." 

In the present analysis, the Green’s function solution" for 
the two-flux differential equations is extended for spectral 
property variations in a composite with general boundary con- 
ditions. This solution for the local radiative heat source is com- 
bined with an implicit finite difference solution of the transient 
energy equation using a tridiagonal matrix. Results are ob- 
tained to validate the numerical method and to illustrate some 
aspects of transient behavior. They compare, for example, gray 
and spectral calculations, and show how additional scattering 
in one layer of a composite can alter the temperature distri- 
bution in an adjacent layer. The solution can be applied to 
examine many transient effects in a two-layer composite with 
spectral properties. 



Fig. 1 Geometry and nomenclature for transient radiation and 
conduction in a two-layer semitransparent composite with spec- 
trally dependent properties and isotropic scattering. 


state is reached, corresponding to the imposed external con- 
ditions. 


Energy Equation and Two-Flux Relations 
for a Two-Layer Composite 

The dimensionless energy equation for transient tempera- 
tures in each layer is 3 


dt(X, t) d 2 t(X, t) _ 1 dq r (X, t) 
dr 1 dX 2 4 6X 

dt(X, t) N 2 d 2 l(X, t) 1 dq r (X, t) 

dr c 21 dX 2 4c 2i dX 


0 <X<5 (la) 

5<X< 1 (lb) 


The last term on the right in each equation is the local heat 
source, including all frequencies, produced by radiative ab- 
sorption, emission, and scattering. Properties are assumed in- 
dependent of temperature. The dq r (X , t )/ dX is obtained in 
terms of the instantaneous temperature distribution by using 
the two-flux relations; these have provided accurate results in 
previous steady-state 13 ' 14 and transient 3 " analyses. With n, and 
n 2 assumed independent of frequency, the spectral two-flux 
equations using the Milne -Eddington approximation are 12 


dq r (X, r, v) _ Ku 
dX 8 

- G,(X, r, v)\ 
dq r (X, r, v) k 2v 


dX 1-8 

- G 2 (X, t, v)] 8<X 


(1 - a iv )[4n 2 e b (X, r, v) 

0 < X< 8 

(1 - n 2 „)[4nle b (X, t, v) 

1 


(2a) 


(2b) 


where G { (X, r, v) and G 2 (X, r, v) are related to q r (X, r, v) 
by 


Analysis 

The analysis is for a plane two-layer composite. Fig. 1, that 
has different spectral radiative absorption and scattering prop- 
erties in each layer. The layers are semitransparent absorbing 
dielectrics with refractive indices larger than one. There is heat 
conduction and isotropic scattering in each layer. The external 
boundaries and the interface between the layers are assumed 
diffuse. The diffuse surface reflectivities are determined from 
the rij using the Fresnel equations; since the tij are assumed 
independent of v, the reflectivities do not depend on frequency. 
Initially, the composite has a temperature distribution T(x, 0), 
and Ti is a characteristic temperature of T(x, 0) used for di- 
mensionless quantities. To begin a transient, the composite is 
subjected to surroundings so that each boundary can receive 
or lose radiant energy, and is in contact with a convecting gas. 
Transient temperature distributions are obtained until steady 


dG t (X, t, v) „ K tv _ 

— =-3-^(X,r,,) 


dG 2 (X, T, V) k 2v _ 

= -3 ~ q r (X, r, v) 


dX 


1 - 5 


0 < X < 5 (3a) 

8 < X < 1 (3b) 


The q r and G are related to the positive ^and negative radiative 
fluxes in Fig. 1 by q r = q7 — q7 and G = 2 (q7 + q r )- 

Boundary, Interface, and Initial Conditions 

For a semitransparent material, radiation is absorbed within 
the material interior and not at a surface, since the surface has 
no thickness. Hence, the relations between conduction and 
convection at the boundaries do not include radiative terms 


dt(X, t) 


dX 


x=o 


Hi 
4 N l 


lA-i - <0. t)1 


(4a) 
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dt(X, t) 
dX 


x=i 


h 2 

4N 2 


Wl, T) - t g2 ] 


(4b) 


At the internal interface, X = 8, there is continuity of temper- 
ature and heat conduction 


g 2 (X, v) 

_ f C l2 (£ 2 , v)cosh e 2v X 4- C 22 (£ 2 , i^sinh e 2v X 8 ^ X < £ 2 
- I C 32 (£>, F)cosh e 2v X 4- C 4 -,(£>, v)sinh e 2v X £> < X ^ 1 

(8b) 


t(8— , t) = ?(84-, t) 


(4c) 


Nt 


dt{X, t) 


ax 


= M 


3r(X, t) 


dX 


(4d) 


The boundary conditions for radiation must include the ef- 
fects of external and internal reflections at the outer surfaces 
that are exposed to spectral fluxes q rl (y) and q r2 { v). The con- 
ditions have been developed by considering the incident, re- 
flected, and transmitted fluxes at a boundary 13 


28 dG x (X, r, v ) 
G,(0, r, v) = 4R,q r ,(v) + - — R 2 - 

3Kj j, oa 


(5a) 


2(1-8) 3G.(X, r, v) 

G 2 ( 1 , t, v) = 4R 5 q r2 ( v) - — R, 


3 k 2 . 


ax 


(5b) 


At the interface between the two layers the relation between 
reflected and transmitted fluxes gives 14 


- o 26 3G,(8, r, v) 

G 2 (8, t, v) =/? 3 G,(S, t, v) 4- R . 4 

3ki„ ax 


(6a) 


To obtain the coefficients C xj (^j, v) to C 4j (^, v), the g,(X, 
v) in Eq. (8) must satisfy the following relations obtained from 
the homogeneous forms of Eqs. (5a) and (5b): 


gi(0, v) - 


26 3g,(X, v) 

3 k Xv 2 ax 


0 


g 2 ( 1, *0 + 


2(1 - 6) 3g 2 (X, is) 

3 k 2 „ 6 ax 


= 0 

X— 1 


(9a) 

(9b) 


By using Eqs. (8a) and (8b) in Eqs. (9a) and (9b) to evaluate 
C„ (£, v) to C 4J ($, v), the g,(X, $, v) reduce to 


Si(X, f„ v) 

f C 21 (£, v)(sinh e,„X 4- P 2v e Xv cosh e^X) 0 ^ X < 

~ 1 C 4 |(£„ v)(sinh e Xv X — F Xv cosh e Xv X) g x < X < 6 

(10a) 


g 2 (X, v) 


CM, 

C 42 (&, 


v)(sinh <? 2 „X — 
r*)(sinh e 2v X — 


/s,, cosh e 2i ,X) 
F 3w cosh e 2v X ) 


5 < X < & 
£ < X < 1 
(10b) 


There is a discontinuity in G at an interface, so that G 2 (5, r, 
v) =£ G,(8, t, v). The radiative flux, however, is continuous at 
the internal interface, and so from Eq. (3) 

_8_ 3G,(X, t, v) 
k Xv ax 

Green’s Function to Determine G(X, t, v) 
in a Two-Layer Composite 

To find dq r (X, r)/3X for Eq. (1), the G,(X, r, v) are needed 
for use in Eq. (2), which is then integrated over all v. Green’s 
functions are developed to obtain Gj (X, r, v) in terms of e b (X, 
r, v). Equations (3a) and (3b) are differentiated in X, and 
3<? r (X, t, v)tdX is then eliminated from Eqs. (2a) and (2b). The 
resulting differential equation for G, (X, r, F) is 


where F Xv = (/? 3 sinh e Xr 8 + P 4v e Xv cosh e Xv 8)/(R 3 cosh e Xy 8 4- 
P Av e Xv sinh e Xv 8), F 2v = cosh e 2 „5/sinh e 2v 8, and F 3v = (sinh 
e 2v 4- Pt, v e 2v cosh <? 2l .)/(cosh e 2v 4- P 6v e 2v sinh e 2v ). 

The C 2 /i ; , *•>) and C 4 ,-(& v ) in Eq. (10) are obtained by 
applying the two conditions 15 that g y (X, §, v) is continuous at 
X = £„ and that for increasing X, the 3g, (X, r-)/aX passes 
through a discontinuity of — 1 at each X = This yields the 
relations for G 21 (^,, v) and C 4 ) (£i, v): 

C 21 (fi, v)(sinh 4- P 2v e Xv cosh e lt ,^) 

= C 4 i(£i, I^)(sinh cosh <?i„£,) (11a) 

C 41 (fi, ^)e lv [cosh e Xv ^ x — F Xv sinh e,^) 

- C->,(£„ v)ei^(cosh 4- P 2v e Xv sinh e x „g x ) = -1 

(lib) 


1 - 6 3G 2 (X. t, v) 


ax 


(6b) 


a 2 G,(X, r, v) 

ax 2 

J = 1. 2 


— e%Gj{X , r, v) = —4ej„n~e h (X, t, r>) 


(7) 


For convenience, the t functional notation is omitted in the 
following derivation. The homogeneous part of Eq. (7) is 
satisfied by cosh(e^X) and sinh(e^X), and so the Green’s func- 
tions have the form 15 


and for C 22 (£>, v ) and C 42 (£ 2 , v): 


C 22 (g 2, v)(sinh e 2l ,£ 2 ~ F 2l/ cosh e 2v g 2 ) 

~ C 42 (^ 2 , ^(sinh e 2v % 2 — F 3u cosh e 2v £ 2 ) (12a) 

C 42 (£ 2 , F)e 2 ^(cosh e 2 „^ 2 — F 3 „ sinh e 2 „i; 2 ) 

— C 22 (€ 2 , i>)e 2v (cosh 2 — F 2v sinh e 2 „g 2 ) = —1 

(12b) 


g,(X, v) 

Cndu i^)cosh e x „X 4- C 21 (^!, r>)sinh e lv X 0 ^ X < ^ 
v)cosh e^X 4- C 4 i(^ l5 r-)sinh e Xv X g x < X S 8 

(8a) 


Equations (11a) and (11b) are solved for C 2 i(£i, v) and C 4i (^i, 
v), and Eqs. (12a) and (12b) for C 22 (§>, v) and C 42 (^ 2 , v). After 
substituting into Eqs. (10a) and (10b), the final Green’s func- 
tions are 


(tanh e Xv £ x + P 2v e x ^) cosh e Xv ^ x 


g.(X, ^ v) = 


Denom lv 

{Ri tanh[e, v (8 - gQ] 4- F > 4y g ly }cosh[g 1 ^(8 - g,)3 
Denom lr 


{i? 3 tanh[^„(5 - Xx)] + P^e ly }cosh[e x ^8 - X)] 0 < X < & 


(13a) 


[(tanh e x „X + P 2 „e ijcosh e Xv X\ % x < X ^ 8 
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£,(X v) = \ 


cosh[e 2l ,(g 2 ~ 5)] 
Denom->„ 


{tanh[e 2v (l — X)] + P 6v e 2v )cosh[e 2l ,(\ — X)] 8 < X < g 2 

cosh|> 2 „(X ~ 6)] £> < X < 1 


{tanh[e 2v (l - fe)] + P 6 t ,g 2 v }cosh[g 2 Xl ~ &)] 
Denom 2 „ 


(13b) 


where 


Denomi„ = e,„{[(/? 3 + / > 2l ,P 4 „e^)tanh e Xv 8 
+ e Xv (P 4v + P 2v /? 3 )]cosh e lt/ 8} 


G 2 (X, t, v) = G /l2 (X, r, v) 

+ 4e 2 „n 2 f g 2 (X, r, g 2 , v)e b (£ 2 , r, *>) dg 2 (15b) 


Denom 2 „ = e 2t ,{ 1 + e 2 „P 6 „ tanh[e 2 „(l — 5)]} 
X cosh[e 2 ,,(l — 5)] 


The gj (X, gj, v) in Eq. ( 1 3) are used to account for the non- 
homogeneous term containing e b (X, r, v ) in Eq. (7) when com- 
puting Gj (X, r, v) at each time_ during the transient. To ob- 
tain the complete solution for G,(X, t, ig), the homogeneous 
solution is needed for Eq. (7). This has the form G hj (X, r>) = 
Zs ljr> sinh <? ;i ,X + E 2jv cosh e jv X. The boundary conditions in 
Eqs. (5) and (6) are applied and the result is four simultaneous 
equations for the E ljp and E 2jl , coefficients: 


e Xv 8 

Kl„ 


P 2t ,e lv E Uv - E 2lv = —4R x q rX (v) 

'e x „8 


cosh c,„S 1 E Uv 4- 


sinh e Xv 8 I E 2iv 


- cosh e3rS ] E„. 

- r ^ 1 ~ g) sinh e,„sl = ~C lv 

L K ** J 


(14a) 


(14b) 


(P 3 sinh e x „8 + P 4v e Xv cosh e Xv 8)E llv 
4- (R 3 cosh e x „8 4- P 4v e Xv sinh e Xv 8)E 2lv 
— (sinh e 2v S)E l2v — (cosh e 2v S)E 22v = C 2 „ (14c) 

(sinh e 2v + P 6v e 2v cosh e 2v )E x2 „ 

4- (cosh e 2 „ 4- P 6y e 2 „ sinh e 2v )E 22v = 4 R 5 q r2 (v) (14d) 


where by using the Green’s functions 


C u, — 4 R 3 e\ v n x 

. 8 f S (tanh + P 2 „<?,„)co sh e,^ x ^ 

X — — e t (i, 

k J 0 Denom tl/ 


v) dg, 


C 2 „ = 4e\ v n\ 

{tanh[<? 2 „(l — 


r*l 

« 6 


X e b (^ 2 , v) dg 2 


gz)] + P 6 „c 2 „}cosh[c 2t ,(l 
Denom 2 „ 


&)] 


The forms for two spectral bands can be written in the same 
way as in Ref. 3, where a spectral Green’s function analysis 
is given for a single layer. 

Numerical Solution 

Starting with a specified f(X, 0), the blackbody function e b (X, 
0, v) is known. The Green’s functions are used in Eqs. (15a) 
and (15b) to obtain Gi(X, 0, v) and G 2 (X, 0, v), and the d 
q r (X, 0, v)/dX is then evaluated from Eq. (2). This spectral 
quantity is integrated over all v, so that it can be used in the 
energy equation that involves total energy quantities (the re- 
sults given later are for two spectral bands with the same cutoff 
frequency in both layers). Using Eq. (1), i(X, t) is integrated 
forward one time increment using the finite difference relations 
that follow. Then the new temperature distribution is used with 
the Green’s functions to obtain G/X, At, v) and dq r (X, At, v)f 
dX. The dq r (X, At, v)!dX is integrated over all v and is used 
in Eq. (1) to continue on to obtain t(X, t) at the next time step. 

From Eqs. (la) and (lb), the following implicit finite dif- 
ference algorithms are used at the interior points of a uniform 
grid to advance the temperature distributions each At in each 
layer: 


N x At 


(AX) 
MA 


2 t(X - AX, t + At) + 


1 + 


27V) At 
(AX) 2 


t(X, t 4- At) 


TKW kx + a * t+ 



f(X - AX, t + At) + 


2A7 2 At 1 

1 + c 21 ( AX) 2 J /( ’ 


(16a) 
t + At) 
At dq r \ 


, t(X + AX, t + At) = f(X, t) 
c 21 ( AX ) 2 V 4c 21 ax 


8 < X < 1 


(16b) 


For the boundary points X = 0 and 1, Eq. (16) is applied, 
and the r(- AX, t + At) and ?(1 + AX, t + At) that appear 
are eliminated by using the boundary conditions Eqs. (4a) and 
(4b), expressed as 

t(AX, t + At) — f(— AX, t + At) 

2AX 


By adding G /y (X, v) to the nonhomogeneous solutions obtained 
by using g,(X, g, v), the general spectral solution of Eq. (7) in 
each layer at each time is 


G,(X, t, v) = G/,,(X, t, v) 


r 

" 4N X Usl 


— f( 0, t + At)] 


t(l + AX, t + At) — ?(1 — AX, t + At) 
2AX 


4- 4 e~ Xv n\ 


*8 

JO 


g x (X, T, g l5 v)c 6 (g), T, V) dg, 


(15a) 


4N 2 


[r(l, r + 


At) - t >;2 ] 
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This yields the finite difference relations at X = 0 and 1 as 
1(0, r + At) 


2NAr HAA 

AX ) 2 2AXJ 


2NAt Ardq r 

(AX? r(AX, T + At) = /(0, T) - — ^ 


H, At 
+ 2AX tgl 

(17a) 


2N 2 At 

c 2 ,(AX) 2 


r(l - AX, r+ At) + 1 + 


2N 2 At HAt 


c 2 i(AX) 2 2c 21 AX 


At dg.. 

X i( 1 , t + At) = t( 1 , t) 


HAt 

+ —t g2 


4 c 2l ax ,. T 2 c 2i AX 


(17b) 


The condition at the internal interface, Eq. (4d), is written 
by using three-point derivatives, and Eq. (4c) is substituted to 
give 


A, 

2 AX 


[1(5 - 2 AX, t) - 4r(5 - AX, t) + 3r(5, t)] 


n 2 

2AX 


[-1(5 + 2 AX, t) + 4r(5 + AX, t) - 3 1 ( 8 , 


r)] 


Equation (16) is written for t at X + AX and X — AX, and 
these relations are used to eliminate the nontridiagonal ele- 
ments t (8 — 2AX, t) and t (8 + 2AX, t). The result is a tri- 
diagonal finite difference expression at the internal interface 


+ ^ 


t( 5 - AX, T + At) 4 2(A, + N 2 )t{ 8 , t + At) 


2i (AX) 2 


- 2 N 2 


t (8 + AX, t + At) 


(AX) 2 


At 

(AX) 2 / dq, 

ax 


[r(S - AX, T) + c 21 i(5 + AX, t)] 


dqr 

ax 


(18) 


The t(X, t) is advanced to the next At using Eqs. (16-18) 
solved as a tridiagonal matrix using the algorithm in Ref. 16. 
As described in Ref. 4, overall energy balances were checked 
at each time step by using the calculated temperature distri- 
butions to evaluate the amounts contributed by radiation, con- 
vection, and internal energy storage. For the grid sizes used, 
the overall energy balance was satisfied within 0.5% through- 
out the transient calculations. After checking various grid 
sizes, it was found that 81 evenly spaced points across each 
layer gave accurate grid-independent solutions for k Jv < 3 in 
each layer. For k Jv > 3, 101 points were used, which gave 
accurate solutions for k Jv £ 10. A time increment of At = 
0.0025 provided stable and time-accurate solutions for the re- 
sults calculated here, as validated by comparisons with some 
previous numerical solutions using the exact radiative transfer 
relations. 


Results and Discussion 

The objective of the analysis was to develop a convenient 
and accurate method for obtaining transient temperature dis- 
tributions in a two-layer semitransparent composite with spec- 
tral radiation properties, isotropic scattering, and differing re- 
fractive indices in each layer. The composite is subjected to 
external convection and radiation on both sides. A Green’s 
function solution was derived for the spectral two-flux radia- 
tive transfer equations in two layers, and for the results given 
here, the solution was applied in two spectral bands with the 


same cutoff frequency in both layers. By summing over the 
two bands, the Green’s function is used to obtain the radiative 
energy source term in the transient energy equation at each 
time step. The energy equation was solved with an implicit 
finite difference method that includes convective boundary 
conditions, and the required temperature and heat conduction 
conditions at the internal interface. Some of the transient re- 
sults that will be presented are to verify the method for spectral 
calculations. This is done by comparing transient temperatures 
with previous spectral steady-state and transient results ob- 
tained or verified with the exact radiative transfer equations. 
The other figures that follow illustrate the variety of transient 
results with spectral effects that can be obtained with the pres- 
ent combined analytical and computational procedure. 

The transient temperatures given here begin from a uniform 
temperature distribution t(X, 0) = 1, although the solution ap- 
plies for a general distribution t(X, 0). For a partial validation 
of the present method, the solid lines in Fig. 2 are transient 
results for a two-band symmetric case using the same param- 
eters as the transient solution in Fig. 6c of Ref. 3. The two- 
band spectral optical thicknesses of the layers in order of in- 
creasing frequency are k 1v = k 2v = 5 and 0.05. For t > 0, the 
composite is exposed on both sides to a high-temperature gas 
at T,i = T g 2 = 4T„ and there is large convective heating on 
both sides. The surroundings remain at T s] = T s2 — T-„ so that 
as the temperatures rise, radiative cooling occurs from within 
the composite interior. Radiative cooling has a strong effect, 
and the transient interior temperatures are always considerably 
below the surface temperatures. At steady state, the interior 
temperatures in the central region of the composite are close 
to the uniform temperature for a composite that is opaque. The 
present results duplicated those in Ref. 3, where a single spec- 
tral layer was analyzed, and additional calculations were then 
made to illustrate the effect of heat capacity in the second layer 
of a two-layer composite. The dashed lines in Fig. 2 illustrate 
how increased heat capacity can delay the transient response. 
Heat capacity, however, does not alter the steady-state tem- 
perature distribution. 

Figure 3 provides a comparison with a numerical solution 
of the exact transfer equations for steady state given in Fig. 3 
of Ref. 14. The transient results at large time agree very well 
with this steady-state solution, which is for a two-layer com- 
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Dimensionless coordinate, X = x/D 


Fig. 2 Verification of two-layer analysis with scattering, and 
transient effect of increased heat capacity in the second layer; 
parameters: q rX = q r2 = 1.0 4 , H 1 = H 2 = 5, t gl = t e2 - 4, 5 = 0.5, and 

c 2 i — 1 . 
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Fig. 3 Comparison of two-layer, two-band transient results at 
large time, with steady-state exact numerical results for radiative 
heating on one side; parameters: fl lt , = fl 2v - 0, q rl = 1.2 4 , q r2 = 
0.25 4 , Hi = H 2 = 1, t sl = 1, t s2 = 0.25, 5 = 0.5, and c 21 = 1. 

posite with each layer having two spectral bands. For these 
results, the composite is heated by blackbody radiation at X = 
0 with temperature T si = 1.27). As the surface temperature rises 
at X = 0, convective cooling occurs, since the gas temperature 
remains at the initial temperature T gl = T, There is stronger 
cooling by both radiation and convection at X = 1, where there 
are lower blackbody surroundings and gas temperatures, T s2 = 
T g2 = 0.257,. For comparison, the steady-state solution is 
shown (dot- dashed line) for a composite that is opaque, and 
so radiation exchange occurs only at the external surfaces, and 
there is only heat conduction within the interior. The transient 
solution has large internal radiation transfer and, at steady 
state, is considerably different from the opaque result. 

After the solution method was verified in Figs. 2 and 3, and 
by some additional comparisons, it was used to examine a few 
transient effects to illustrate how the solution can be applied; 
these results are shown in Figs. 4-7. A few steady-state so- 
lutions are also shown for an opaque composite, and it is ev- 
ident that in these examples, considerable radiation effects are 
present. 

In Fig. 4, the effect of adding scattering to the first layer is 
illustrated as a possible method to protect the second layer 
from external radiant heating. For the solid lines, the two-band 
spectral optical thicknesses of the first layer in order of in- 
creasing frequency are K iv = 2 and 0.1. For the dashed lines, 
a uniform amount of scattering has been added to both bands 
so k ]v - 6.9 and 5. The first layer has thus been made more 
opaque, and it tends to shield the second layer from the radi- 
ation incident at X = 0; this significantly reduces the temper- 
atures in the second layer. During the transient, however, the 
maximum temperature is increased in the first layer. At steady 
state, the dashed results approach somewhat more the temper- 
ature distribution for an opaque composite. 

To illustrate the importance of spectral effects, the results 
without scattering in Fig. 4 are compared in Fig. 5 with tem- 
peratures when each layer is gray. In each layer the larger of 
the two-band optical thicknesses in Fig. 4 was used, and so, 
in Fig. 5, k,„ = 2 and k 2v = 5 for all v. The bands have been 
eliminated that are more transparent, so that radiation incident 



Fig. 4 Influence of scattering in first layer on transient temper- 
atures after exposure to radiative heating on one side; parame- 
ters: q,\ — 2 4 , q r 2 = l 4 , H x ~ H 2 — 4, t gI = t g2 — 1, S = 0,5, c 2t — 1, 
and KjT s i = 4000 ftm K; solid lines, f = 0 with = 2 
and 0.1; dashed lines, H,„ = 0.7101 and il 2 „ = 0.98 with = 6.9 
and 5. 



Fig. 5 Comparison of two-band spectral results with those for a 
gray composite; parameters: = 0, q ri = 2 , q r2 ~ 1 * H t 

= H 2 = 4, t gl = t g2 = 1, 8 = 0.5, and c 21 = 1; solid lines, k 1v = 2 and 
0.1, k 2v = 5 and 0.5, and \ CJ T sl = 4000 /xm K; dashed lines, #c,„ = 
2 and k 2v - 5. 

at X = 0 is not transmitted as readily into the composite. In- 
cident radiation is absorbed strongly in the region at small X. 
The result is increased maximum temperatures that are closer 
to X - 0, and reduced temperatures in the second layer. By 
comparison with the dot— dashed line in Fig. 4 it is evident 
that the optical thicknesses in Fig. 5 are significantly smaller 
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Dimensionless position, X = x/D 


Fig. 6 Effect of thermal conductivity in the first layer on tran- 
sient temperature distributions; parameters: = il 2v = 0, q rl = 

2 4 , q r2 = l 4 , Hi = H z = 4, t s i = t g2 = 1, S = 0.5, c 21 = 1, and A cj T sl = 
4000 /nm K. 



Fig. 7 Transient temperatures after exposure to radiation on one 
side without convective cooling; parameters: (I,, = fl 2v = 0, q rl = 
2 4 , tf r2 = l 4 , t gI = t g2 = 1, 5 = 0.5, c 21 = 1, and KjT sl = 4000 K. 


than those required for the composite to behave like an opaque 
material. 

The results in Fig. 6 illustrate the effect of heat conduction 
in the first layer, which is the layer adjacent to the radiant heat 
addition. Except for the Nj values, the parameters are the same 
as for the solid lines in Fig. 4. For the solid lines, = N 2 = 
0.1, as compared with A/j = AL = 0.015 in Fig. 4. The increased 
Nj have a substantial effect in reducing the temperature gra- 
dients near the external boundaries, and the transient distri- 
butions are more uniform as expected. The dashed and dot- 


dashed curves show the effect of increasing and decreasing the 
conductivity of the first layer. An increase to N t = 1 removes 
practically all of the temperature peak near X — 0. Decreasing 
the conductivity to N, = 0.01 produces a substantial tempera- 
ture peak in the first layer and a large gradient near X - 0 
throughout the transient. 

As a final illustrative example. Fig. 7 examines the effect of 
removing the convective cooling at the X — 0 boundary. The 
solid lines are the same as the solid lines in Fig. 6, where H\ 
- 4; for the dashed lines, H, - 0. As expected, removing con- 
vective cooling at X = 0 raises the temperatures near this 
boundary throughout the transient, and at X = 0 the tempera- 
ture distributions have a zero derivative. At steady state, the 
entire temperature profile is raised because of the reduced 
cooling on one side. The straight lines for an opaque composite 
at steady state with a uniform thermal conductivity illustrate, 
as in Fig. 4, that there are substantial internal radiation effects. 

Conclusions 

A combined analytical and numerical method was developed 
to obtain transient temperatures in a composite of two semi- 
transparent layers with spectral radiation properties and with 
isotropic scattering. A Green’s function solution was derived 
for the spectral two-flux equations in the composite. Summing 
over all frequencies provides the local radiative heat source. 
This was combined with an implicit numerical solution of the 
transient energy equation. The analysis accounts for heat con- 
duction, and internal radiant absorption, emission, and isotro- 
pic scattering. An advantage of the method is that isotropic 
scattering is included without additional complexity. Starting 
from an initial temperature distribution, the layer is heated or 
cooled by external radiation and convection. The resulting 
steady-state and transient temperature distributions were veri- 
fied by comparing with results where the radiative energy 
source was numerically evaluated from the exact spectral ra- 
diative transfer equations. Illustrative examples demonstrate 
how the method can be used to examine various heat transfer 
effects, and they show the comparison of transient spectral and 
gray calculations. An example is shown where the addition of 
scattering to a layer is used to shield a second layer from 
sudden exposure to radiative heating. The effects of changing 
the heat conductivity in part of the composite and removing 
the convective cooling at one boundary are also illustrated. 
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